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INTRODUCTION 


The  purpose  of  this  document  is  to  present  the  mathematical  treat¬ 
ments  used  in  the  description  of  the  various  aspects  of  a  computer  simu¬ 
lation  of  a  Space  Based  Radar  Surveillance  (SBRS)  system,  bearing  in  mind 
that  the  validity  of  the  results  of  a  simulation  rests  on  the  conditions  of 
validity  of  the  mathematical  models  describing  the  various  phenomena 
involved.  We  will  therefore  make  a  special  effort  to  establish  as  clearly 
as  possible  the  conditions  of  validity  of  the  equations  and  the  limitation 
of  the  various  models  that  we  used. 

Our  objective  in  making  a  computer  simulation  of  a  SBRS  system  is  to 
be  able  to  simulate,  as  accurately  as  possible,  the  interaction  of  a  SBRS 
system  with  various  attack  scenarios.  The  duration  of  these  attacks  is 
limited  to  a  few  hours.  From  that  consideration  comes  our  first  assumption: 
we  will  neglect  the  movement  of  the  Earth  around  the  Sun.  The  Earth  will  be 
considered  as  having  a  fixed  center  of  rotation. 

MOVEMENT  OF  THE  SATELLITES 


The  motion  of  the  satellite  is  assumed  to  obey  the  unperturbed 
Kepler  laws.  It  means  that  the  Earth  is  represented  by  a  homogenous 
sphere  with  its  mass  concentrated  at  the  center.  The  orbit  of  the  satel¬ 
lites  are  then  perfect  ellipses.  No  perturbations  of  any  kind  are  included 
in  the  calculations  of  the  satellite  position.  The  true  anomaly  f  is  then  a 
function  of  E,  the  eccentric  anomaly  (1): 


f  3  sin-l(/-.L-e.2  sjn.E)  (1) 

1-e  cos  £ 

But  the  eccentric  anomaly  E  is  a  transcendental  equation  of  time: 

E  -  e  sin  E  *  (t  -  T)  (2) 

l 3 

where:  T  =  integration  constant  determined  by  specifying  the  time  origin. 

ii  =  k2m 

k2  =  gravitational  constant 
m  =  mass  of  the  Earth 
i  -  semilatus  rectum 
e  =  eccentricity 


We  solve  the  equation  2  by  using  the  Brown’s  method  as  described  in 
reference  1.  The  true  anomaly  was  then  calculated  from  equation  1.  The 
radial  distance  r  is  given  by 


r 


1 

1  +  cos  f 


(3) 


The  movements  of  the  satellite  in  the  orbital  plane  is  then  fully  described 
by  the  polar  coordinates  (r,  f)  as  Illustrated  in  Figure  1. 


We  set  the  time  origin  such  that  for  T  =0  the  satellite  is  at  its 
perigee.  We  also  set  a  system  of  rectangular  coordinate  ERE  whose  origin  is 
at  the  center  of  the  Earth.  The  X  axis  is  pointing  toward  the  perigee  and 
the  Z  axis  is  perpendicular  to  the  orbital  plane.  The  satellite  turns  in 
the  positive  direction.  The  direction  of  rotation  of  the  satellite  may  be 
inverted  by  inverting  the  direction  of  the  Z  axis  (see  Figure  2).  The 
orientation  of  the  X,Y  and  Z  axis  of  the  ERF  system  is  independant  of  time 
since  the  rotation  of  the  Earth  around  the  Sun  is  neglected. 

It  is  also  convenient  to  define  two  other  systems  of  coordinats 
having  their  origin  at  the  center  of  the  Earth.  The  first  system  is  called 
the  ERF  system.  Its  Z  axis  is  pointing  toward  the  North  Pole  and  its  X  axis 
is  in  the  equatorial  plane  of  the  Earth  at  a  longitude  of  0*  at  and  only  at 
time  T  =  0.  The  orientation  of  the  X,Y  and  Z  axis  of  the  ERF  system  is  not 
a  function  of  time  (see  Fig.  3). 

The  second  system,  the  ERR  system,  is  similar  to  the  ERF  system 
except  that  the  X'  axis  is  always  pointing  at  longitude  0,  in  such  a  way 
that  the  system  is  rotating  with  the  angular  speed  of  the  Earth.  The  XYZ 
and  X'Y'Z'  systems  are  illustrated  in  Figure  3.  The  geographical 
coordinates  (EGR)  of  latitude  and  longitude  are  illustrated  in  Figure  4. 

A  third  system,  the  ORR  system  is  used  to  display  the  results  on  a 
Tektronix  screen.  It  is  a  system  that  is  topocentric  to  some  point  of 
observation  defined  by  the  operator.  That  point  of  observation  is  rotating 
with  the  Earth.  The  origin  of  the  ORR  system  is  at  the  point  of 
observation.  The  X  axis  is  pointing  away  from  the  center  of  the  Earth,  the 
Z  axis  is  pointing  northward  and  the  Y  axis  is  always  in  a  plane  parallel  to 
the  equatorial  plane  (see  Fig  5).  Because  of  the  definition  of  the  system 
of  coordinates  the  Y  and  Z  components  only  are  displayed  on  the  screen. 

In  order  to  find  the  position  of  the  satellite  at  a  certain  time  in 
the  ORR  system,  the  position  of  the  satellite  is  first  calculated  in  the  ERE 
systems,  then  those  coordinates  are  transfered  to  the  ERF  system,  and 
finally  to  the  ORR  system,  taking  into  account  the  Earth  rotation  when 
necessary.  The  transfers  are  implemented  using  a  matrix  of  rotation 
calculated  by  the  routine  ANGROT  and  proper  translations. 

A  set  of  parameters  that  describes  the  orbit  of  the  satellite  without 
redundancy  or  ambiguity  in  a  way  that  is  independant  of  time  has  been 
selected.  Those  parameters  are  stored  in  an  array  called  SYSAT  whose  format 
is  described  in  Figure  6.  The  direction  cosine  are  calculated  in  the  ERF 
system.  The  user  of  the  simulation  can  change  the  content  of  SYSAT  at 
will. 

COVERAGE  AREAS 

We  build  in  our  simulation  the  two  types  of  coverage  areas  associated 
with  the  electronic  scanning  mode  (see  Fig.  7)  and  with  the  mechanical 
scanning  mode  (see  Fig.  8).  The  blind  spot  centered  at  nadir  has  been  set 
at  a  grazing  angle  of  60*.  The  outer  limit  of  the  area  is  defined  in  the 
eleventh  column  of  SYSAT.  The  width  of  the  mechanically  scanned  area  is 
either  set  in  SYSAT  (workspace  CURAR)  or  calculated  from  the  radar 
parameters  (workspace  SCAN). 
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FIGURE  3  -  RELATIONSHIP  BETWEEN  THE  FIXED  AND  ROTATING  SYSTEM  OF 
COORDINATES  (ERF  AND  ERR  SYSTEMS) 


FIGURE  4  -  GEOGRAPHICAL  SYSTEM  OF  COORDINATES  (EGR) 
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The  size  of  the  footprint  on  the  ground  is  defined  by  the  angular  width  of 
the  beam  and  the  geometry  of  the  system  as  illustrated  in  Figure  16.  The 
relationship  between  the  parameters  are  given  by  the  equations  4,5,6,  7,8 
and  9  applied  respectively  to  the  far  and  the  near  border  of  the 
footprint.  The  user  has  to  choose  the  amount  of  overlap  of  successive  foot 
prints.  A  few  examples  are  shown  on  Figure  17.  The  distance  between  the 
center  of  the  footprints  is  given  by 


0 


/ 100 
50 


=-S)  x  R. 


where  0  is  the  amount  of  overlap  in 


Of 


•O  m 


We  have  assumed  here  that  two  successive  footprints  have  the  same  size  and 
that  the  curvature  of  the  Earth  is  negligible  when  we  consider  two 
successive  footprints. 


MOVEMENT  OF  AIRPLANES 


The  movement  of  the  airplanes  are  defined  and  controlled  by  the  data 
stored  in  the  array  FP  whose  format  is  illustrated  in  Figure  18.  To  every 
airplane  corresponds  one  plane  of  information  in  FP.  The  first  line  of  each 
plane  contains  the  initial  parameters  of  the  airplane: 

-  altitude,  longitude  and  latitude  of  the  initial  position  of  the 
ai rplane 

-  the  Initial  horizontal  speed,  heading  and  climb  rate  of  the 
airplane. 

-  the  current  time  that  is  reset  to  0  before  every  run  of  the 
simulation 

-  the  start  time  of  the  airplane. 

The  second  and  third  line  contains  the  data  necessary  to  define  a 
change  of  flight  plan  at  a  certain  time  breakpoint.  They  contain 

-  the  new  horizontal  speed,  heading  and  climb  rate, 

-  the  time  of  the  breakpoint 

The  rest  of  the  line  does  not  have  to  be  predefined  by  the  operator 
it  is  used  to  store  the  data  during  the  calculations  (see  Figure  19). 

The  fourth  line  contains  the  termination  time  of  the  flight  plane 
(that  is  defined  by  the  operator)  and  will  contain,  after  the  execution  of 
the  program,  the  final  position  of  the  airplane. 

When  the  program  is  run,  it  reads  the  flight  plan  stored  in  FP, 
calculates  the  position  of  the  airplane  at  the  required  time  accordingly  to 
the  flight  plan  and  stores  the  current  time  and  the  position  of  the  airplane 
in  FP.  If  the  current  time  is  between  the  start  time  and  the  time  of 
breakpoint  1,  this  is  done  in  line  two  of  the  appropriate  plane;  if  the 
current  time  is  between  the  time  of  breakpoint  1  and  the  time  of  breakpoint 
2,  the  storage  is  done  in  line  3,  etc...  So,  at  the  end  of  the  execution, 
the  final  position  of  the  airplane  is  stored  in  line  4  (see  Figure  20). 
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FIGURE  15  -  FORMAT  OF  THE  RADPAR  ARRAY 


FIGURE  14  -  ANGULAR  DEFINITION  OF  A  BARRIER  (BSR  SYSTEM) 


FIGURE  12  -  FORMAT  OF  THE  ARRAY  BARRIER 


-  15  - 


parameters  contained  in  RADPAR  and  SYSAT.  More  details  on  the  scanning 
process  may  be  found  in  the  next  section. 


The  parameters  that  define  the  location  of  the  barrier  are  stored 
in  an  array  called  BARRIER  (see  Figure  12).  A  barrier  is  set  up  by 
the  location  of  the  two  ends  of  its  spine,  the  arc  of  great  circle  AB  (see 
Figure  13).  The  width  defined  in  the  array  BARRIER  is  the  width  at  both 
ends  of  the  barrier.  The  X3  axis  should  be  going  through  the  center  of  the 
cord  AB.  If  we  now  consider  a  spherical  system  of  coordinate  (BSR)  with  the 
angle  P  measured  from  the  Z3  axis  (see  figure  14),  the  barrier  will  be 
def i ned  by 


Lg  <  0  <  Lg 


L*  <  (f  -  *)  <  L< 


(10) 


where  L0  and  L.  are  0  and  $  half-angles  subtended  by  the  barrier  from  the 
origin  of  the  system  of  coordinates.  Then,  a  program  calculates  the 
location  of  the  four  corners  of  the  barrier.  The  boundaries  of  the  barrier 
are  the  areas  of  great  circle  joining  the  four  corners.  To  check  if  a  scan 
area  or  a  coverage  area  has  some  overlap  with  a  barrier,  one  has  only  to 
express  the  coordinates  of  the  borders  of  that  area  in  the  BSR  system  and  to 
check  if  some  of  them  are  within  the  limits  of  the  barrier. 

It  is  convenient  to  define  two  systems  of  coordinates  associated  with 
the  barrier.  These  systems  will  be  used  to  determine  if  some  part  of  the 
scan  lines  or  of  the  coverage  areas  of  the  radars  are  within  the  barrier. 

The  first  system  (BRR)  is  a  rectangular  system  centered  on  the  Earth.  The  Z 
axis  (Zo  on  Figure  13)  Is  parallel  to  the  line  AB  joining  the  extremities  A 
and  B  of  the  arc  AB  that  is  the  spine  of  the  barrier.  The  X  axis  (Xn  on 
Figure  13)  passes  through  the  middle  of  the  line  AB.  The  second  (BSR) 
system  (see  Fig.  14)  Is  a  spherical  system  derived  from  the  BRR  system  with 
the  angle  $  measured  from  the  Z  axis. 

SCANNING  PROCESS 


The  size  of  the  scan  line  generated  by  a  radar  depends  on  the 
parameters  of  the  radar.  They  are  stored  in  an  array  called  RADPAR  whose 
format  is  illustrated  on  Figure  15.  The  width  of  the  beam  is  defined  as  the 
50%  point.  The  angular  width  of  the  beam  is  defined  from  the  antenna  size 
through  the  equation: 

0  ,  Lll.  (il.) 

11  x  S 


where 


\  3  wavelength  of  the  radar 
s  3  size  of  the  antenna. 
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FIGURE  9  -  GEOMETRICAL  PARAMETERS  INVOLVED  IN  THE  CALCULATION  OF 
THE  COVERAGE  AREA 


The  geometry  of  the  situation  is  illustrated  in  Figure  9.  AID  is 
the  grazing  angle  stored  in  SYSAT.  A  few  useful  relationships  between 
the  parameters  involved  are  listed  here.  There  is  no  approximation  in 
these  equations  and  the  Earth  curvature  has  been  taken  into  account.  The 
Earth  is  assumed  to  be  spherical  with  a  radius  R. 

0  =  sin*1  (5—  x  sin(-2  +  ALD)) 
HS  2 

(4) 

BRC  =  HS  x  tg  0 

(5) 

MRC  =  tg© 

(6) 

vnr  _  (-BRC  x  MRC)  +  /  R2  +  (MRC 

x  R)2  -  BRC2 

(7) 

1  +  MRC2 

RPC  =  (MRC  x  XRC)  +  BRC 

(8) 

RDMAX  =  /  (HS  -  XRC ) 2  +  RRC2  (9) 

The  angle  0  can  be  either  the  angle  associated  with  the  outer  border 
of  the  coverage  area  or  the  angle  associated  with  the  inner  border  of  the 
coverage  area.  The  parameter  RDMAX  is  used  in  the  program  that  checks  if 
there  is  an  airplane  within  the  coverage  area.  RDMAX  is  the  maximum 
distance  from  the  satellite  where  you  may  find  an  airplane  within  the 
coverage  area  (on  the  same  side  of  the  Earth  as  the  satellite). 

SYSTEM  OF  BARRIERS 

Out  of  many  SBRS  strategies,  we  choose  to  simulate  the  one  in 
which  the  radars  on  the  satellites  operate  only  when  the  scan  area  overlaps 
some  predefined  areas  called  "barriers".  Most  of  the  time,  in  our 
simulation,  the  radars  are  in  "search  mode".  It  means  that  the  scan  areas 
of  the  radars  are  calculated  at  fixed  intervals  and  a  check  is  made  to  know 
if  there  is  an  overlap  with  the  barrier.  If  a  new  overlap  is  found,  the 
program  will  run  a  routine  to  find  the  time  of  the  first  contact  with  the 
barrier  and  from  then,  the  radar  will  operate  in  scan  mode  as  long  as  its 
scan  area  is  not  out  of  the  barrier.  The  scan  area  of  the  electronic  and 
mechanically  scanned  radar  are  shown  on  Figures  10  and  11  respectively.  The 
width,  the  length  and  the  location  of  the  scan  areas  are  calculated  from  the 
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FIGURE  7  -  ELECTRONIC  COVERAGE 


FIGURE  6  -  FORMAT  OF  THE  ARRAY  SYSAT 
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5  -  SYSTEM  OF  COORDINATE  TOPOCENTRIC  TO  THE  POINT  OF 
OBSERVATION  (ORR) 
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FOOT  PRINT 


FIGURE  16  -  GEOMETRY  OF  THE  FAR  AND  NEAR  LIMITS  OF  THE  FOOTPRINT 


100%  OVERLAP 


FIGURE  17  -  ILLUSTRATION  OF  PERCENTAGE  OF  OVERLAP 


FIGURE  18  -  FORMAT  OF  THE  FP  ARRAY 
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FIGURE  19  -  REPARTITION  OF  THE  DATA  IN  FP  BEFORE  AN  EXECUTION 
OF  A  FLIGHT  PLAN 


FIGURE  20  -  REPARTITION  OF  THE  DATA  IN  FP  AFTER  AN  EXECUTION 
OF  A  FLIGHT  PLAN 
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This  particular  approach  to  the  handling  of  the  data  related  to 
the  movement  of  the  airplanes  has  first  been  proposed  and  implemented  by 
Or.  8.  Bridgewater.  We  adapted  the  concept  to  the  particular  needs  of  our 
simulation. 

It  is  to  be  noted  that  the  system  can  operate  only  for  positive 
time  progression.  The  program  will  give  erroneous  data  if  the  time  is  run  in 
reverse. 

The  position  of  the  airplane  is  always  calculated  from  the  preceeding 
breakpoint  or  from  the  start  point.  The  airplanes  are  assumed  to  move  on 
great-circles.  The  new  position  of  the  airplane  is  first  calculated  in  a 
system  of  coordinate  (PRC)  centered  on  the  Earth  with  the  Z  axis 
perpendicular  to  the  plane  containing  the  great-circle  in  which  the  airplane 
is  moving.  The  X  axis  passes  through  the  initial  position  of  the  airplane 
(see  Figure  21).  The  angular  movement  of  the  airplane  and  its  new  position 
are  calculated.  The  position  of  the  airplane  is  then  calculated  in  a  system 
centered  on  the  Earth  with  a  topocentric  orientation  of  its  axis  by  using  a 
proper  matrix  of  rotation  MRX  (calculated  by  the  program  VALMFP).  The 
position  of  the  airplane  is  then  calculated  for  display,  in  the  system  ORR 
by  using  the  matrix  of  rotation  MRFO. 

RELATIVE  SPEEDS 


The  first  step  to  establish  the  signal  to  noise  ratio  associated  with 
a  particular  satellite-target-Earth  geometry  consists  of  calculating  the 
relative  speed  of  the  satellite  and  the  target  and  the  relative  speed  of  the 
satellite  and  the  clutter  from  the  Earth  behind  the  target  as  seen  by  the 
satellite.  We  also  build  a  routine  that  calculates  the  relative  speed 
between  the  satellite  and  II  sample  points  located  on  the  scan  line  (see 
Fig.  22). 

We  first  calculate  the  speed  of  the  satellite  and  of  the  target  in 
the  ERF  system  by  using  the  matrix  of  rotation  already  defined.  The 
component  of  the  vector  joining  the  satellite  and  the  target  or  the 
satellite  and  the  sampling  point  on  the  scan  line  in  the  ERF  system  to  get 
the  line  on  which  the  speed  will  be  projected.  The  satellite,  the  target 
and  the  Earth  speed  are  then  projected  on  the  proper  lines  and  the  Doppler 
shifts  are  calculated. 

ROTATION  OF  COORDINATES 


In  order  to  perform  the  many  rotations  of  coordinates  generated  by  the 
programs,  we  build  a  routine  that  generates  the  matrix  of  rotation  needed  to 
go  from  a  X'Y'Z'  system  to  a  XYZ  system  and  the  reverse.  The  routine  ANGROT 
works  from  the  components  of  the.  X'  and  V  vector  in  the  XYZ  system  and 
calculates  the  matrix  MR  to  go  from  the  X'Y'Z'  system  to  the  XYZ  and  the 
inverse  matrix  IMR  to  go  from  the  XYZ  to  the  X'Y'Z'  system  (see  Figure  23). 
One  may  be  surprised  by  the  complexity  of  ANGROT  that  contains  more  than 
4000  different  cases.  One  should  realize  by  example,  that  each  time  for, 
X'Y'or  Z'are  in  one  of  the  planes  XOY,  XOZ  and  XOY,  we  get  a  degenerated 
situation  that  has  to  be  handled  in  its  own  way.  A  close  examination  of  the 
problem  leads  to  the  identification  of  a  very  large  number  of  such 
situations.  This  is  why  a  check  procedure  has  been  added  to  ANGROT.  It 
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FIGURE  21  -  SYSTEM  OF  COORDINATES  RELATED  TO  THE  AIRPLANE  MOVEMENT 
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FIGURE  22  -  ILLUSTRATION  OF  THE  RELATIVE  SPEED 


FIGURE  24  -  MAP  PRODUCED  ON  A  TEKTRONIX 
SCIENCE  ASSOCIATES  INC.  SOF' 


simply  submits  the  X'  axis,  as  defined  in  the  XYZ  system,  to  an  IMR 
rotation.  If  the  rotation  is  right,  we  will  get  as  an  answer  the  X'  axis  as 
defined  in  the  X'Y'Z'  system  (1  00). 


MAP  PACKAGE 

It  is  also  possible  to  add  to  our  display  a  nice  outline  of  the 
boundaries  of  the  continents,  as  seen  from  a  point  of  observation  (see 
Figure  24).  The  software  package  generating  the  data  is  written  in  Fortran 
language  was  gratuitously  given  to  our  group  by  Decision  Science  Associates 
Inc.  This  package  is  designed  to  calculate  the  projection  of  the  boundaries 
of  the  Earth's  continental  features  for  various  types  of  projection. 

Mr.  B.  Cook  adapted  these  programs  to  our  computer  system  and  generated,  in 
a  proper  format,  a  file  containing  the  data  associated  with  an  orthonormal 
projection  for  an  origin  located  at  90°  of  longitude  West  and  45°  of 
latitude  North.  That  file  was  read  in  APL,  transferred  to  an  APL  workspace 
and  processed. 

We  generated  only  the  data  necessary  to  produce  a  projection  centered 
at  90"  West  and  45*  North.  Other  maps  may  be  generated  at  will  by  running 
the  programs  with  other  sets  of  parameters. 


CONCLUSION 

A  brief  description  of  the  mathematical  treatment  involved  in  SBRS 
simulation  has  been  presented.  The  principle  of  operation  of  the  routine 
was  also  described.  Special  attention  was  paid  to  the  statement  of  the 
condition  of  validity  and  of  the  approximation  present  in  the  various 
models. 
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